﻿using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using UFPE.CIn.Abc.Primality.Util;
using System.Diagnostics;

namespace UFPE.CIn.Abc.Primality.Algorithms
{
    public class LucasSelfridge : BaseAlgorithm
    {
        public override string AlgorithmName
        {
            get
            {
                return "Lucas-Selfridge";
            }
        }


        private static LucasSelfridge instance;

        public static LucasSelfridge Instance
        {
            get
            {
                if (instance == null)
                {
                    instance = new LucasSelfridge();
                }

                return instance;
            }
        }

        //***********************************************************************
        // Implementation of the Lucas Strong Pseudo Prime test.
        //
        // Let n be an odd number with gcd(n,D) = 1, and n - J(D, n) = 2^s * d
        // with d odd and s >= 0.
        //
        // If Ud mod n = 0 or V2^r*d mod n = 0 for some 0 <= r < s, then n
        // is a strong Lucas pseudoprime with parameters (P, Q).  We select
        // P and Q based on Selfridge.
        //
        // Returns True if number is a strong Lucus pseudo prime.
        // Otherwise, returns False indicating that number is composite.
        //***********************************************************************

        /*public bool IsPrime(BigInteger val)
        {
            cron = Stopwatch.StartNew();

            BigInteger thisVal = val.AbsoluteValue();

            switch (thisVal.PreliminaryChecks())
            {
                case BigInteger.Status.None:
                case BigInteger.Status.Composite:
                    cron.Stop();
                    timeLastRun = cron.ElapsedMilliseconds;

                    return false; //goto ReturnFalse;
                case BigInteger.Status.Prime:
                    cron.Stop();
                    timeLastRun = cron.ElapsedMilliseconds;

                    return true; //goto ReturnTrue;
            }

            return Helper(thisVal);
        }*/

        public override bool Helper(BigInteger thisVal)
        {
            // Do the test (selects D based on Selfridge)
            // Let D be the first element of the sequence
            // 5, -7, 9, -11, 13, ... for which J(D,n) = -1
            // Let P = 1, Q = (1-D) / 4

            long D = 5, sign = -1, dCount = 0;
            bool done = false;
            bool isPrime = false;

            while (!done)
            {
                int Jresult = BigInteger.Jacobi(D, thisVal);

                if (Jresult == -1)
                    done = true;    // J(D, this) = 1
                else
                {
                    if (Jresult == 0 && Math.Abs(D) < thisVal)       // divisor found
                        goto Return;

                    if (dCount == 20)
                    {
                        // check for square
                        BigInteger root = thisVal.sqrt();
                        if (root * root == thisVal)
                            goto Return;
                    }

                    //Console.WriteLine(D);
                    D = (Math.Abs(D) + 2) * sign;
                    sign = -sign;
                }
                dCount++;
            }

            long Q = (1 - D) >> 2;

            /*
            Console.WriteLine("D = " + D);
            Console.WriteLine("Q = " + Q);
            Console.WriteLine("(n,D) = " + thisVal.gcd(D));
            Console.WriteLine("(n,Q) = " + thisVal.gcd(Q));
            Console.WriteLine("J(D|n) = " + BigInteger.Jacobi(D, thisVal));
            */

            BigInteger p_add1 = thisVal + 1;
            int s = 0;

            for (int index = 0; index < p_add1.dataLength; index++)
            {
                uint mask = 0x01;

                for (int i = 0; i < 32; i++)
                {
                    if ((p_add1.GetDataAt(index) & mask) != 0)
                    {
                        index = p_add1.dataLength;      // to break the outer loop
                        break;
                    }
                    mask <<= 1;
                    s++;
                }
            }

            BigInteger t = p_add1 >> s;

            // calculate constant = b^(2k) / m
            // for Barrett Reduction
            BigInteger constant = new BigInteger();

            int nLen = thisVal.dataLength << 1;
            constant.SetDataAt(nLen, 0x00000001);
            constant.dataLength = nLen + 1;

            constant = constant / thisVal;

            BigInteger[] lucas = BigInteger.LucasSequenceHelper(1, Q, t, thisVal, constant, 0);

            if ((lucas[0].dataLength == 1 && lucas[0].IsZerothBitEqualToZero()) ||
               (lucas[1].dataLength == 1 && lucas[1].IsZerothBitEqualToZero()))
            {
                // u(t) = 0 or V(t) = 0
                isPrime = true;
            }

            for (int i = 1; i < s; i++)
            {
                if (!isPrime)
                {
                    // doubling of index
                    lucas[1] = thisVal.BarrettReduction(lucas[1] * lucas[1], thisVal, constant);
                    lucas[1] = (lucas[1] - (lucas[2] << 1)) % thisVal;

                    //lucas[1] = ((lucas[1] * lucas[1]) - (lucas[2] << 1)) % thisVal;

                    if ((lucas[1].dataLength == 1 && lucas[1].IsZerothBitEqualToZero()))
                        isPrime = true;
                }

                lucas[2] = thisVal.BarrettReduction(lucas[2] * lucas[2], thisVal, constant);     //Q^k
            }


            if (isPrime)     // additional checks for composite numbers
            {
                // If n is prime and gcd(n, Q) == 1, then
                // Q^((n+1)/2) = Q * Q^((n-1)/2) is congruent to (Q * J(Q, n)) mod n

                BigInteger g = thisVal.gcd(Q);
                if (g.dataLength == 1 && !g.IsZerothBitDiffOfOne())         // gcd(this, Q) == 1
                {
                    if ((lucas[2].GetLastData() & 0x80000000) != 0)
                        lucas[2] += thisVal;

                    BigInteger temp = (Q * BigInteger.Jacobi(Q, thisVal)) % thisVal;
                    if ((temp.GetLastData() & 0x80000000) != 0)
                        temp += thisVal;

                    if (lucas[2] != temp)
                        isPrime = false;
                }
            }

        Return:
            return isPrime;
        }
    }
}
